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^ ■ We study the old problem of the uniqueness of chemical evolution models by ana- 

. lyzing a set of multiphase models for the galaxy NGC 4303 computed for a variety of 

(—] ! plausible physical input parameters. 

We showed in Molla et al. (1999) that multiphase chemical evolution models for 
O ■ the three Virgo cluster galaxies NGC 4303, NGC 4321 and NGC 4535, were able to 

' reproduce the observed radial distributions of spectral indices Mg2 ad Fe52. Chemical 

■ evolution models may, however, fit the present-epoch radial distributions with differ- 

ent star formation histories, thus we need to include time-dependent constraints to 
^ . the problem. The two spectral indices above depend on the time-averaged history of 

star formation, but are in turn affected by the well known age-metallicity degeneracy 
which prevents the disentangling of age and metallicity for stellar populations, another 
uniqueness problem we also discuss. 

Our aim is to determine if the possible input parameters may be strongly constrained 
when both radial distributions (nebular and stellar) are used. In order to accomplish 
this we run a large number of models (500) for NGC 4303 varying the model input 
parameters. Less than 4 % of the models (19) fit the present day observational data 
within a region of 95% probability. The number of models reduces to a ~ 1 % (6) 
when we also ask them to reproduce the time-averaged abundances represented by the 
spectral indices. Thus, by proving that only a small fraction of the models are favored 
in reproducing the present day radial abundance distributions and the spectral indices 
data simultaneously, we show that these spectral indices provide strong time-dependent 
additional constraints to the possible star formation and chemical histories of spiral 
disks. 
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Subject headings: galaxies: abundances - galaxies: evolution- galaxies: spiral - 
galaxies: stellar content 



1. Introduction 

Variations of abundances with galactocentric radius have been widely observed in external 
spiral galaxies (Diaz 1989; Skillman et al. 1996; Garnett et al. 1997), mostly via observations of Hll 
regions which represent the abundances in the gas phase at the present time. Since Tinsley (1980), 
a large number of chemical evolution models have been developed, by including different physical 
mechanisms (see Gotz & Koppen 1992; Koppen 1994) in order to reproduce these radial gradients. 

In this work we use a multiphase chemical evolution model (MCEM) where an infall dependent 
on galactocentric distance and a two-step star formation law, whereby molecular clouds are formed 
before stars, are the basic hypotheses. The outcome is star formation simulating a power law in 
the gas density with an exponent n > 1. This model was first applied to the Solar Neighborhood 
(Ferrini et al. 1992), to the whole Galactic (MWG) disk (Ferrini et al. 1994, hereafter FMPD), 
and to the bulge (Molla & Ferrini 1995). Next the same model was applied to disks (Molla et al. 
1996), and bulges (Molla et al. 2000) of a set of spiral galaxies of different morphological types, the 
observed radial distributions of diffuse and molecular gas, oxygen abundances and star formation 
rate being successfully reproduced. It seems therefore that the MCEM produces results in good 
agreement with the data representing the present-day interstellar medium. 

It is well known, from studies of the MWG, however, that by tuning in an appropriate way both 
the history of star formation and the infall rate it is possible to reproduce the same present-time 
radial distributions using very different evolutionary scenarios (see, for examples and further refer- 
ences, Molla et al. 1992; Tosi 1996). Tosi (1988) early addressed this problem of the uniqueness of 
theoretical chemical evolution models. In fact, in most cases differences in this kind of models have 
no, or very small, effect on the observed radial distributions as was shown in Tosi (1996), where a 
comparison of models built by different groups was discussed. The only important characteristic 
that discriminates among the existing chemical evolution scenarios is the resulting time evolution 
of the radial gradients of abundances. This time evolution depends on the ratio SFR/Infall, which 
becomes the crucial quantity in determining whether abundance gradients steepen or flatten with 
time. If infall is the source of star-forming gas, the gradient flattens as the galaxy evolves and 
consumes this gas. Conversely, it steepens with the time if the effect of the infall gas is mostly 
that of diluting the metallicity. In other words, present time Hll abundances are not sufficient by 
themselves to discriminate among the possible models of disk evolution. Since chemical abundance 
data other than those derived from H II regions are almost inexistent for external galaxies, a check 
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on the time evolution of chemical models is only possible, in principle, for the MWG, where ob- 
servational data for individual stars of different ages exist. This comparison was performed for the 
MCEM in Molla et al. (1997) 2. 

An independent test for the evolutionary history of external spiral disks as predicted by the 
chemical evolution models was proposed by Molla et al. (1999, -hereafter Paper II). There we made 
use of the radial distributions of spectral indices Mg2 and Fe52, obtained by Beauchamp &; Hardy 
(1997), hereafter Paper I, and Beauchamp (1997) for three galaxies in the Virgo Cluster, NGC 4303, 
NGC 4321 and NGC 4535. The MCEM applied to these galaxies, reproduced the present radial 
distributions and also predicted, for the first time in spiral disks, the main characteristics of the 
radial distributions of Mg2 and Fe52 indices. The reader must remember, as pointed out above, 
that the calibrated spectral indices represent abundance time-averages over the galaxy chemical 
history, as opposed to the accumulated effect shown by the extreme Pop. I objects such as Hii 
regions or the B stars. 

But even in the case under discussion we face an uniqueness problem because the spectral 
indices synthesis models also show a degeneracy in the plane Mg2-Fe52 (or (Fe)) plane: Indeed, age 
produces almost the same effect as metallicity thus compromising their simultaneous determination 
for a given single stellar population, if only these two quantities are used (Barbuy 1994; Worthey 
1994; Buzzoni et al. 1992, 1994). Other spectral indices, mostly those related to the Balmer 
lines, must be used in order to maximize the age-metallicity discrimination and to constraint the 
evolutionary parameters (Freitas-Pacheco 1998; Delisle 1998). The latter studies refer to elliptical 
galaxies where a high degree of coevality is assumed. Spiral disks are subject to a further difficulty 
in that stellar disks are composite systems and the time vagaries of star formation must be taken 
into account. They are on the other hand richer in information and physical constraints since data 
on the gas abundance and density are available and must be reproduced too. 

Thus, there exists a possible ambiguity in that different star formation histories in disks are 
possible that would lead to the same final averaged spectral indices and/or present gas distributions. 
Are there alternatives to the SFR histories found in Paper II which may reproduce at the same 
time the present gas radial distributions and the stellar content information? Do both problems 
- uniqueness and degeneracy - go in the same or in opposite directions? In other words, can we 
constrain the possible sets of evolutionary histories on the basis of all the available observations? 

In this work, we attack this problem by modifying the input parameters of the MCEM for 
the galaxy NGC 4303 for which we compute a large number (500) of different chemical evolution 
models. Then, by using a optimization, we choose those able to fit the present-time nebular 
observations with probabilities larger than 97.5 %. Next, radial distributions of Mg2 and Fe5270 are 



^Another possible way to test the time evolution of the models is of course to calculate the evolution of abundances 
with redshift and to compare the predictions with the observations. This was successfully accomplished for damped 
Lyman alpha systems (Ferrini et al. 1997), but this procedure does not permit examining the time evolution of radial 
gradients, or the abundance variations in different regions of the same disk. 
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computed by using the chemical abundances and star formation histories predicted by the above- 
selected evolutionary models as input for the synthesis procedure employed to calculate spectral 
indices. We will determine if the comparison of these spectral indices distributions with the data 
helps in the final selection of the evolutionary histories valid over the disk. 

In Section 2 we briefly describe the MCEM. The computed models as well as the method of 
selection among them, based on their predictions for the present epoch are also in Section 2. The 
calculation of synthesis models and the radial distributions of the spectral indices Mg2 and Fe5270, 
as well as their comparison with observations, are shown in Section 3. We discuss these results in 
Section 4. Finally, our conclusions are presented in Section 5. 

2. The Multiphase Chemical Evolution Model (MCEM) 

2.1. A Brief Description 

We remind the reader that the appropriate check of validity of the basic multiphase model, 
here adopted, has already been done in previous work, leaving open however the issue of uniqueness. 
Therefore, we only summarize here the MCEM whose detailed properties, as applied to spiral disks, 
are described in Ferrini et al. (1992, 1994); Molla et al. (1996) and in Paper II of this series, and 
concentrate on the comparison among models. 

In the MCEM a protogalaxy is assumed to be a spheroid composed of primordial gas with 
total mass M(R). For NGC 4303 the latter is calculated from the rotation curves obtained either 
via radio (Guhathakurta et al. 1988) or optical (Distefano et al. 1990) observations. The galaxy 
is divided into concentric cylindrical regions 1 kpc wide. The model calculates the time evolution 
of the halo and disk components belonging to each cylindrical region. The halo gas falls into the 
galactic plane to form the disk, which in the multiphase framework is a secondary structure formed 
by the gravitational accumulation of the diffuse gas, ^ij, at a rate /. The infall rate / is inversely 
proportional to the collapse time scale TcoU, which is assumed to be dependent on galactocentric 
radius through an exponential function with a scale length Xd: 

TcoiiiR) = Toe-(^-^°)/^° (1) 

The characteristic time scale tq above is defined for every galaxy as being that of a region 
located at a radius Rq, equivalent to the Solar Neighborhood (hereafter SN) in the MWG, which is 

itself scaled from Rcff = 62.6 arcscc (Henry (1992)). Thus, the adopted radius Rq for this galaxy is 
6 kpc. The value for the characteristic time scale may be determined, using the total mass of the 
galaxy, through the expression: 



To oc Mg ^/^T 



(2) 
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from Gallagher et al. (1984), where Mg is the total mass of the galaxy in 10^ Mq and T is its 
age, assumed 13 Gyr. Thus, we obtain tq for each galaxy, via the ratio: 

ro = ro(M9,5aVM9,MWG)-V2 (3) 

where Mg^MWG is the total mass of MWG and the value of r© = 4 Gyr, corresponding to the 
Solar Neighborhood, was determined in Ferrini et al. (1992). 

In the various regions of the galaxy (halo or disk) we allow for the different phases of matter 
aggregation: diffuse gas (g), clouds (c, only in the disk), low-mass (si,m < 4:Mq) and massive stars 
(s2, m > 4:Mq), and remnants. The mass in the different phases changes through several conversion 
process, all of them related to the rate of the corresponding star formation process: 

1. Star formation by the gas spontaneous fragmentation in the halo (oc Kg]^^) 

2. Cloud formation in the disk from diffuse gas (oc ng^'^) 

3. Star formation in the disk from cloud-cloud collisions (oc He?) 

4. Induced star formation via massive star-cloud interactions (oc acs2) 

5. Diffuse gas restitution from these processes 

The rates for these processes are proportional to the parameters K, /x, H, and a, which, in turn, 
depend on galactocentric radius through the equations derived in FMPD. The proportionality 
factors in these equations are the corresponding efficiencies of the various processes, that is the 
efficiency of the halo star formation, ek, and, in the disk, the probability of cloud formation, e^^, of 
cloud — cloud collision, eu, and of the interaction of massive stars with clouds, €«. 

The term is assumed constant for all halos, thus being independent of morphological type. 
The last term in the list above, containing the induced star formation, is associated to local processes 
and, as a result, its coefficient is considered independent of both position and morphological 
type. However, the other two rates, and eh, depend on the Hubble type and/or the arm class. 
Therefore they must be chosen for each galaxy, their variation range being determined following 
the arguments given by Ferrini & Galli (1988), as discussed in Molla et al. (1996) and in Paper II. 
These efficiencies resulted larger for earlier morphological types and lower for the later ones, and 
this is useful for the selection of initial values for a given galaxy as we showed in Paper II. Within 
this range, a fine tuning selection for these efficiencies and en is performed by choosing those for 
which the present time radial distributions (those of oxygen abundances, star formation rate, and 
the atomic and molecular gas surface densities) are well reproduced. We must stress that although 
each galaxy has its own set of efficiencies, the corresponding parameters ji, H, and K vary along 
the galactocentric radius. 
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The adopted initial mass function (IMF) is taken from Ferrini et al. (1990). The enriched 
material is the result of the restitution of processed material by dying stars and depends on their 
nucleosynthesis processes, their IMF (and delayed restitution) and their final fate, either via quiet 
evolution or via Type I or II supernova (SN) explosions. Nucleosynthesis yields used here are from 
Renzini k, Voli (1981), for low and intermediate mass stars, and from Woosley & Weaver (1995) for 
massive stars, respectively. The type I supernova explosions release mostly iron, following Nomoto 
et al. (1984) and Brach et al. (1986), at a slower rate than type II, with the result that the iron 
appears at least 1 Gyr later than the a-elements (Oxygen, Magnesium ...) ejected by the massive 
stars. 

The initial conditions are such that when star formation begins early in the halo zone and 
later in the disk (a secondary structure in our model), most of the barionic galactic mass has 
already accumulated. Boundary conditions are arbitrarily chosen: the first and simplest choice is 
to consider galaxies as closed systems. This choice does not strongly constraints the internally- 
driven phenomena - SN winds for example are not as prevalent during spiral galaxies evolution 
as to produce significant out-gassing, although they can influence internal redistributions. The 
environment in which galaxies live may have a sensibly larger importance (Henry 1993, 1996; 
Pilyugin & Ferrini 1998). We discuss here only closed models, keeping in mind that this is a first 
order approximation. 

In the same spirit, we do not include at present in our models the effect of bars. It is quite 
evident, even from elementary dynamical considerations, that barred galaxies must present flatter 
radial gradients, as it has indeed been observed (Edmunds &; Roy 1993; Martin &; Roy 1992, 1994, 
1995). We will simply consider that an eventual bar will influence the radial distributions only in 
the few internal kpc. 

Summarizing, we will retain the same basic scenario applied to other spirals models but will 
vary the model input parameters, namely: the characteristics collapse time scale, tq, according to 
the total mass, and the two efficiencies, and en, depending on the morphological type. The 
combination of selected parameters will therefore be different for each galaxy. 

2.2. The Computed Models 

Following the arguments above we had estimated in Paper II a value ^ 8 Gyr as the char- 
acteristics collapse time scale and a scale length Ad of 4 kpc, for the galaxy NGC 4303. After 
estimating these parameters defining the collapse law, we chose the efficiencies e^^ and en, (0.22, 
0.01), as those able to reproduce the observed radial distributions for the atomic gas (Hi, Warmels 
1988; Cayatte et al. 1990), the molecular gas (H2, Kenney Sz Young 1988, 1989), the star formation 
rate surface densities obtained from Ha fluxes from Kennicutt (1989), and flnally the gas oxygen 
abundances (Henry 1992; McCall et al. 1985; Shields et al. 1991). In the present work will use the 
oxygen abundances as were corrected by Pilyugin (2000), which seem more adequate. 
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We proceed now to change the input parameters in order to take into account the uncertainties 
in the selection of that first set. Thus, for instance, the total mass depends on the maximum rotation 
velocity which, in some cases is different when calculated either from radio Hi or optical data. For 
example, for NGC 4303, Vmax = 150 km/s increases to a value around 216 km/s, from the first to 
the second data set. This factor of ~ 1.5 in velocity would translate into a factor of ~ 2.25 in tq. 
Therefore, we must increase or decrease the collapse time scale by running models with 5 possible 
values for the characteristics collapse time scale: r = 1, 4, 8, 12 and 16 Gyr. 

The relation between collapse time scale and mass makes it evident that TcoU must change 
with galactocentric radius. If the total mass surface density follows an exponential form such as 
the surface brightness of spiral disks does, then the resulting TcoU must increase along the radius 
with a scale length Ad oc Re (where Re is the scale- length for the surface brightness distribution). 
This Re varies with the wavelength band and, moreover, it usually decreases for later types of 
galaxies and is larger for the earlier ones. Therefore, although we have so far used in our standard 
models Ad = 4 kpc for all spiral disks, we must keep in mind the possibility of variations of this 
parameter. Thus, we take other 5 possible values for this scale length: Ad = 1, 4, 8, 12 and 16 kpc. 

The effects of changing efficiencies or the collapse time scale are not equivalent: by increasing 
the cloud and star formation efficiencies it is easier to form molecular clouds and also to destroy 
them to produce stars, thus resulting in a star formation rate with higher absolute values, but 
maintaining the same star formation history. In other words, the maximum of the star formation 
rate will occur always at the same time for models with the same TcoU- Alternatively, a longer 
collapse time scale means a lower infall rate on the disk, resulting into a star formation rate which 
shifts its maximum to a later time, producing therefore a slower evolution even with the same 
absolute value of its maximum. Prom Eq.l the scale length Ad, changes the rate of evolution 
among different radial regions, increasing the differences among them if it is small or decreasing 
them when it is large, thus acting upon the radial gradients of abundances, which become steeper 
(or flatter) when the parameter Ad is shorter (or larger). 

In order to check the effect of the degeneracy problem in the plane Mg2 — Fe52, it is crucial 
to analyze the differences among the ages of stellar populations predicted by models with varying 
collapse times or length scales, and their effects on spectral indices. 

Efficiency values and collapse times are constrained variables that change together. These 
probabilistic efficiencies may take values from up 1, and we assume that they depend on the 
morphological type, following the equations: 



,-T2/15 



(4) 



and 



en = e-^V5 



(5) 
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where T is the Hubble type (see Molla et al. 2002), which gives values similar to those from 
Paper II where T = 4, the morphological type assigned to NGC 4303. 

Because these values might change if we vary the infall rate law, we have also run models 

with 20 possible values for each set (e^,,e//), calculated from the above expressions, simulating 20 
different morphological types from to 10 (by including intermediate values as 0.5, 1.5, etc.). 
Thus, we ran a total of 500 models, with all possible variations of these three parameters, tq, Xd 
and the set (e^,€i^) -or equivalently the value of T-. 



We must now select which set of parameters among these 500 models are adequate to represent 
the observed present-time radial distributions of abundances, star formation rate, and diffuse and 
molecular gas density. 

In order to do that, and following the technique used by Tosi (1988) in her uniqueness problem 
study, we compute the indicator, which measures the proximity of a model to the region limited 
by observational error bars: 



where ^caiji is the computed quantity, in each radial region i, l^b^^i is the corresponding observed 
quantity at the same galactocentric distance, and (j, is the error in quantity. 

We compute our models for radius from 2 up to 14 kpc in steps of 1 kpc, by avoiding the 
inner region or bulge because it has a different geometry and structure. Therefore, i takes values 
from 2 to 14. The value o"j is taken as ~ 0.2 dex for the oxygen abundances, the current error in 
the estimates of abundances, and ~ 2Mq/pc^ for the Hi surface density, a mean value estimated 
for these types of data. For the star formation rate and the molecular cloud surface density we 
assume errors as larger as 40% and 50%, respectively, due to the larger uncertainties involved in 
the derivation via Ha and CO fluxes. In particular, we must take into account that the H2 masses 
in a given region arc estimated through a scaling factor a which depends on the metallicity of this 
region (Verter & Hodge 1995; Wilson 1995). Due to the uncertainties involved in the calculation 
of this parameter, the corresponding radial distribution of (t(H2) are considered doubtful (Taylor 
k Klein 2001). 

After computing the indicator for each radial distribution of oxygen, star formation and gas 
surface densities produced by our 500 models we must choose the possible models or the range in 
the input parameter which may produce models with results falling within the confidence region. 
In order to do that, we first select the model corresponding to the lowest values of x^, Xmin' ^'^^ 
each set of observations, that is the maximum-probability model for each constraint. 



2.3. The most probable models from the present day observations 




(6) 
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We now select those models for which — Xmin < A, A being the selected value for a given 
probability which takes into account the free parameters of the sample and selects the confidence 
region for the parameter set From the number of radial regions, parameters and observational 
constraints, we search for the parameter region with — Xmin < 5.62, that is probabilities larger 
than 97.5% for each individual data set or larger than 90% for all of them over the region of 
superposition. In this way, for each model and observational constraint, we have a value of x^, and 
its corresponding probability P. 

In order to choose the best models, we must obtain the region of parameters that minimize 
their x^ for each radial distribution of data separately, and then find the region of coincidence. 
Accordingly, we represent these regions of high probabilities for each observational constraint, all 
on the same figure so as to delimit the values of the most probable input parameters. When there 
are two free parameters in the theoretical model finding the plane parameterl-parameter2 for which 
each set of constraints has a large probability of being well fitted, is relatively easy. In practice 
we are dealing with a 3 parameters space-ro, Xd and the set (e^ and en), or its equivalent T 
value. Thus, in order to find the region of parameters that reproduce all observations with a high 
probability we will relay on projections of these planes taking two parameters at a time. 

In fig. 1, (panels a to e), we show the probability contours for each radial distribution plotted, 
as function of type T and of the collapse time scale tq (in Gyr). These 5 panels, one for each Ad 
value, are the successive projections on the plane T — tq of our 3-dimensional regions along the 
Ad axis. Those parameters which give probabilities larger than 97.5% for the oxygen abundance 
12 + log{0/H) (PoH > 0.975) are represented by solid lines. The dotted lines represent the star 
formation rate surface density (Psfr > 0.975) and the atomic and molecular surface densities 
( Phi > 0.975 and > 0.975), are represented by the long and short dashed lines. The 
shaded regions are the coincidence regions where all constraints are reproduced simultaneously 
with P > 90%. 

In panel a) there are no successful models, all of them produce probabilities smaller than 97.5% 
of fitting the observational data. It is therefore clear than A > 4 kpc. The panel with the largest 
number of acceptable possibilities correspond to A = 8 kpc where the region of coincidence is the 
largest of the five. The longer the collapse time scale the larger must be the scale length. In effect, 
if Ad = 4 kpc, 4 < tq < 12 Gyr; if Ad = 8 kpc tq > 8 Gyr, while if Ad = 12 or 16 kpc, r > 10 - 12 
Gyr. For all these possibilities, efficiencies fall within a very limited region of 4 < T < 5. 

The latter statement is also apparent from panels f) to j), where the same kind of probability 
contours are shown for the plane T-Ad. Prom them, we can eliminate the shortest collapse time 
scale (r =1): only if r > 4 Gyr it is possible to have coincidence at a level larger than 90 % for all 
data. If r = 4 Gyr, only a value of Ad = 4 kpc is possible. If r = 8 Gyr, then 4 > Ad > 8 kpc. 
Only for larger collapse time scales, larger values of Ad, up to 16 kpc, are also possible. 

^From the graphs above, it is evident that only efficiencies corresponding to a T between 4 to 
5 (i.e., only intermediate and late Spirals) are possible. Selection of other efficiency values would 
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Fig. 1. — Probability contours representing regions where P > 97.5%, for each observational 
constraint. Top panels, from a) to e) as a function of T and tq. Central panels, from f) to j), as 
a function of T and A^. Bottom panels, from 1) to o), as a function of tq, and A^. The solid line 
represents the contours for the oxygen abundances. The same for the star formation rate surface 
density are the dotted lines, and for the atomic and molecular gas surface densities are the long- 
dashed and the short-dashed lines, respectively. The shaded regions (or the full dots when there is 
no region) are the zones of coincidence of these 4 contours, showing the range of parameters more 
probable than 90%. 
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produce poor fittings, mostly in the gas radial distributions, -long and dashed lines-, for which 
increases abruptly if T < 4 or T > 5. The bottom panels from k) to o), represent the plane Tcoi — A 
for T values from 3.5 to 5.5 For values out of this range, there are again no coincidence regions. It 
is also clear from these panels that tq > 4 — 8 Gyr, and A^j > 4 kpc are limiting values to obtain 
good models. Furthermore, it is also evident that the efficiencies values are strongly constrained by 
the gas densities, and that the use of both the diffuse and molecular gas densities, instead of the 
total gas mass or the gas fraction, proves superior in finding models fitting the present- day data. 

The allowed morphological types reduce only to 2 possible values, T = 4 or T = 4.5, in excellent 
agreement with the type T = 4 given for NGC 4303. If T = 5 only a model (with r = 16 Gyr and 
Xd = 8 kpc, is able to fit all constraints. Other T values, or different efficiencies of diffuse gas and 
molecular cloud formation have null probabilities of fitting the gas radial distributions. Since the 
morphological type is at least partially defined by the amount of gas on the disk we do not regard 
this as a surprising result, but are reassured that our efficiencies for the different morphological 
types are well chosen. 

To summarize, the present method yields parameters for which the model results are within the 
region defined by the error bars with a confidence level A. We reduce our possibilities to a Ad larger 
than 4 kpc, a characteristic collapse time scale between 4 and 16 Gyr, and efficiencies corresponding 
to morphological types T= 4 or 4.5. Only 19 models out of 500 satisfy these conditions, thus 
reducing the possible chemical evolution models to a ~ 4% of the initial number. 

We show in Table 1 smd the corresponding probability P for each one of the observational 
constraints, for these 19 models which have probabilities larger than 97.5% to be close to the 
observations (which in turn corresponds to 90% for the combination). There we have in column 
(1) the number identifying our model. Column (2) is the value of T. The characteristics collapse 
time scale, tq, and the scale length, Xd, are in columns (3), and (4). Efficiencies and en are 
given in columns (5) and (6). In Columns (7) to (10) the values of x^, corresponding to the radial 
distributions of abundances, star formation rate, and diffuse and molecular gas surface densities, 
are given. Finally, the probabilities of these distributions to be in a region around the minimum 
value of x^ are in columns (11) to (14). We conclude there is only a model with tq = 4 Gyr and 
only a model with T = 5. All the others correspond to T = 4 or 4.5 with tq > 8 Gyr. 

2.4. Present Day Radial Distributions 

Following the procedure described above we have found 19 different models for N4303 satisfying 
the condition that the radial distributions of abundances, star formation rate and gas densities all 
have probabilities larger than 90% of fitting the observational data. 

We represent all models in Fig. 2 where we see that all of them are falling within the regions 
defined by the observational limits, obtained with the data and the error bars and represented by 
the thick solid lines. We regard the latter subset of models as equally probable within the errors 
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Table 1. Parameters of the Selected Models (P > 90%) 



Model 


T 


r 


Ad 


Em 




2 

XOH 


2 

XSFR 


2 

Xhi 


XH2 


PoH 


PSFR 


Phi 


Ph, 


129 


4. 


4. 


4.5 


.259 


.017 


0.957 


6.040 


0.378 


2.318 


1.000 


0.983 


1.000 


1.000 
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Fig. 2. — Predicted Radial distributions for models which satisfies the condition P> 90% simul- 
taneously for all present time data: a) diffuse or atomic gas Hi surface density; b) molecular gas 
H2 surface density; c) Oxygen abundance as 12 + log{0/H) and d) Star Formation Rate surface 
density. Models are represented by different lines and symbols following the figure. Observational 
data are the full dots while the thick lines limit the region within observational error bars used in 
our realizations. 
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and seek now additional constraints to isolate the best models among them, representing the true 
evolutionary history of the galaxy. 

3. Using Evolutionary Synthesis Models 

Unlike Hii abundances, whose abundances reproduce the cumulative effect of all past gener- 
ations of stars, stellar indices contain information on the evolutionary history of a galaxy in the 
form of an average over time and over the SFR. They provide therefore additional constrains to 
our models and we show in this section how to retrieve this information. 

3.1. The method of calculation 

As pointed out, for N4303 (and two other galaxies of paper II) in addition to the typical 
quantities necessary to define the galactic evolution (i.e., gas, metals,...), we have at our disposal 
measurements of radial distributions of the spectral indices Mg2 and Fe52. The fact that the 
models selected in the preceding sections represent well the present-day radial distributions, does 
not imply that the evolutionary history given by the multiphase model are unique, as no time- 
dependent information has been used so far. In fact, we have obtained results in agreement with 
the present gas data either with fast evolution (such as Model 129 with tq = 4 Gyr) thus producing 
stellar populations biased toward older ages and an earlier enrichment, or with a slower creation 
of stars (such as Model 450 tq =16 Gyr) thus producing metals at a more constant rate. Both 
alternatives yield almost the same present radial distributions. Would they also produce similar 
radial distribution of stellar spectral indices? Wc will now compute these indices for each of our 
models in an attempt to answer this question, and will try to check if the radial distribution of 
spectral indices are equally well reproduced by the adopted set of models. 

We proceed in a self-consistent way as follows. Using the chemical evolution models already 
described we calculate the integrated mass of stars created in a given time interval by the gas, 
and the mean abundance reached at that epoch by the same gas. We consider the resulting stellar 
populations residing in every galactocentric region as the superposition of a set of single stellar 
populations or generations, each one of which is defined by its age and its metallicity. We remind 
the reader that our models provide abundances for 15 elements, including Oxygen, Magnesium and 
Iron, so we know in addition the abundance ratio [Mg/Fe] for each generation. We have in principle 
the necessary ingredients to evaluate observable spectral indices via a synthesis model. 

Spectral index features are calculated with the same method described in (Molla & Garcfa- 
Vargas 2000) by using the Padova group isochrones. The technique in question was developed 
to be used specifically in spiral disks or in regions where a combination of stellar populations 
of different ages and/or metallicities is necessary, and it works by simultaneously computing the 
spectral features and the continuum flux. An isochrone is assigned to each of our stellar generations 



-15- 



using the total abundance Z reached at this time step and the corresponding age. Then, an 
index value is calculated for each star belonging to that generation following its stellar gravity 
log (7, effective temperature T^s and metallicity or abundance, through the fitting functions, and 
a synthetic spectra from atmosphere models is selected. By adding all stellar contributions, the 
contribution of every generation to both synthetic single continuum and line fluxes is computed 
at the same time. Next, continuum and line fluxes are added for all generations, weighted by the 
stellar mass created in every time step, according to the star formation history provided by the 
evolution model. This calculation is performed for all disk regions, generating spectral indices along 
the galaxy disk. 

We use the fitting function from Worthey et al. (1994) for assign the index Fe52 to each star. 
However, the fitting function for Mg2 must depend on the abundance [Mg/H], not [Fe/H], contrary 
to what is assumed in that work. The selection of [Mg/H] or [Fe/H] in spiral disks is important due 
the different rate in the ejection of these elements and the differences among the disk regions along 
the radius. Alternatively, Borges et al. (1995) obtained a fitting function for the index Mg2 by 
including a dependence on [Mg/H]. However this dependence is derived from the relative abundance 
[Mg/Fe]. We prefer to use directly the dependence on [Mg/H]. Therefore, by using the same the 
stellar library from these authors, we have derived our own fitting function which gives us the 
direct dependence of Mg2 on [Mg/H] . The abundance of [Mg/H] is used to calculate Mg2 and the 
abundance of [Fe/H] to calculate Fe52. 



3.2. Spectral indices results 

We have computed in this way the radial distributions of spectral indices for all models of 
the galaxy NGC 4303 with larger than 90% probability in all observational constraints for the 
present-time data, as displayed in Table 1. 

Can the comparison with spectral index radial distributions help us further in selecting the 
best evolutionary model or models, or limit the range of parameters of these?. 

We compute once again the statistical indicator but this time for each one of the two 
distributions as given by spectral indices. The results of these calculations are shown in Table 2, 

where we give the goodness of fit for each model. For each model number given in Column (1), 
Column (2) gives for the Mg2 data. Column (3) for Fe52. The corresponding probabilities are 
in Columns (4) and (5). The radial gradients obtained for each model and by the observational 
data distributions are also given in Columns (6), (7) and (8). 

In Fig. 3 we represent the same panels 1) and m) from Fig. 1 but limiting those resulting 
regions with these new constraints. In panel a), T= 4, we see that there is only a point, this one 
corresponding to Model 129, which satisfies our confidence condition of 97.5% for all observational 
constraints. For panel b), T= 4.5, there still exists a little region shaded by both vertical and 
horizontal lines (that one shaded only by vertical lines correspond to the region obtained with 
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Table 2. Results from synthesis models. 



Model 


XMg2 


XFe62 


PMg2 


PFeB2 


Radial Gradients 


Number 










VO/H 


VMg2 


VFo52 


* 129 


5.30 


3.17 


1.00 


1.00 


-0.079 


-0.0062 


-0.094 


* 228 


6.69 


3.94 


1.00 


1.00 


-0.076 


-0.0048 


-0.074 


* 229 


6.60 


4.89 


1.00 


1.00 


-0.086 


-0.0050 


-0.078 


248 


11.06 


3.49 


0.97 


1.00 


-0.050 


-0.0030 


-0.045 


* 249 


10.51 


4.41 


0.98 


1.00 


-0.059 


-0.0037 


-0.056 


* 329 


8.30 


6.75 


1.00 


1.00 


-0.091 


-0.0048 


-0.075 


348 


11.80 


5.02 


0.95 


1.00 


-0.055 


-0.0031 


-0.047 


* 349 


10.70 


6.29 


0.98 


1.00 


-0.063 


-0.0035 


-0.054 


368 


12.81 


4.90 


0.91 


1.00 


-0.047 


-0.0023 


-0.034 


369 


13.41 


4.82 


0.88 


1.00 


-0.054 


-0.0021 


-0.031 


388 


13.41 


4.82 


0.88 


1.00 


-0.043 


-0.0021 


-0.031 


389 


13.17 


6.01 


0.90 


1.00 


-0.050 


-0.0028 


-0.042 


448 


12.19 


6.04 


0.94 


1.00 


-0.058 


-0.0030 


-0.046 


449 


14.38 


10.92 


0.83 


0.90 


-0.066 


-0.0037 


-0.057 


450 


14.38 


10.92 


0.83 


0.90 


-0.072 


-0.0037 


-0.057 


468 


13.13 


5.94 


0.90 


1.00 


-0.050 


-0.0025 


-0.037 


469 


12.49 


7.41 


0.93 


0.99 


-0.057 


-0.0030 


-0.047 


488 


13.61 


5.87 


0.87 


1.00 


-0.046 


-0.0023 


-0.034 


489 


12.73 


7.32 


0.92 


0.99 


-0.053 


-0.0029 


-0.043 


OBS 










-0.078 


-0.0086 


-0.073 
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Fig. 3. — Probability contours representing regions with P> 97.5 % for all observational constraints 
for T= 4 (panel a) and T= 4.5 (panel b) as a function of r and A. Lines represent the same than 
panels 1) and n) from Fig. 1, while filled squares are the models obtained in Section 3.2 
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chemical evolution models) showing the possible input parameters which reproduce all data with 
the same goodness. There are more than one model, but we have restricted the possibilities to a 
very reduced zone. 

The best model fitting the observational data with P = 1 for all constraints , corresponds to 
Model 229, with a morphological type T = 4.5, a collapse time scale TcoU = 8 Gyr and = 4 
kpc, a model similar to Model B presented in Paper II as the best model. Only 5 other models, 
marked with an *, have also able to fit these new constraints with probabilities larger than 97.5%. 
In particular. Model 129, the only one among the 19 with tq = 4 kpc, have also high values of 
probabilities. All the other possible models, mostly in the low part of the table, that is with 
slow collapses and/or low star formation efficiencies, although reproduce the present day radial 
distributions, are far away from our confidence region for one or both spectral indices. Thus, 
Model 450, with P> 0.99 in Table 1 decreases its probabilities, when the comparison of model 
predictions with data is performed through the spectral indices Mg2 and Fe52, to 0.83 and 0.90, 
respectively. Therefore, we reduce to 6 the possible models from the initial 500, this is only 1 % of 
them are retained by the procedure. 

We may answer to the stated question, can the comparison with spectral index radial distribu- 
tions help us to select the appropriate evolutionary model? with an affirmative statement. There 
are some sets of physically plausible input parameters that provide an adequate enough fit to the 
abundance data. However, these models can not reproduce equally well the spectral indices radial 
distributions. 

Wc show in Fig. 4 the predicted radial distributions of Mg2 (in panel a), and Fe52 (in panel 
b), respectively, for these 6 models. In the graphs of Fe52, observational points in the outer disks 
(R > 10 kpc) of all galaxies tend to increase with radius; we think that some of these points at 
large and moderate radii may be affected by large systematic uncertainties in the sky-subtraction 
procedure. The comparison of models with these Fe52 outer data is not significant. 

For both indices there is a inner region where the observed values raise significantly. This 
region seems to contain a stellar bar which cannot be simulated with our present models, and as 
discussed, we have not included these regions in our calculations. However, we have added the 
results obtained for the bulge or central region with the Model B from Paper II in Molla et al. 
(2000). The only difference between this bulge model and the disk models used here resides in the 
geometry and structure of this region. The input parameter for this bulge model are very similar 
to Model 229. Taking into account that the most of differences among models appears in the the 
outer disk predictions, we can assume that by using the 5 other model input parameters in the 
bulge model would give similar results for this bulge or central region. These results are shown 
in Fig. 4 as open squares located in a mean radius for the bulge, assumed as 1 kpc. We have not 
computed the complete model with the full radial range included, for which a model allowing the 
calculation of radial gradients in the bulge is necessary. We believe that these central values fall on 
the right absolute level, and that the shape of the spectral indices distributions in the inner disk 
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Fig. 4. — Spectral indices radial distributions for NGC 4303. Panel a) the index Mg2 in magnitudes. 
Panel b) The index Fe52 in A. Data are from Beauchamp (1997) and Beauchamp & Hardy (1997). 
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may be also well reproduced. 

Summarizing, the uniqueness issue associated the chemical evolution models is not strong, 
because the models were already reduced to a 5% with just the present day data as constraints. But 
when we use in addition the spectrophotometric indices we limit even more the possible evolutionary 
histories. 

3.3. The origin of the radial variations of Magnesium and Iron 

What is the explanation for the strong observed Mg2/Fe52 radial gradient? Indeed, computing 
the Mg2/Fe52 ratio from the data in columns 6,7 and 8 of Table 2 the gradient turns much more 
steeper in the observational data than predicted by the models. We have no straightforward 
explanation for this behavior. 

From a chemical evolution point of view Magnesium appears early in the evolution, as it is 
produced by massive stars which evolve rapidly. Iron, instead, is created by the explosions of SN-I 
which need a time delay to appear. Although abundances tend to a saturation level, saturation 
is reached earlier for Magnesium, and as a result the Iron radial gradient subsists longer than the 
Magnesium one. 

When stellar averaged abundances arc computed, taking into account that star formation rate 
proceeds more rapidly in the inner disk than in the outer regions, an abundance flattening results 
which is very similar for both indices. In fact, the averaged stellar relative abundance < [Mg/Fe] > 
is almost zero for our 6 models, without appreciable variations along the disk. 

We must remember, however, that spectral indices are not directly equivalent to abundances 
because they also depend on the mean age of the parent stellar populations. Therefore, one might 
think that the explanation for the strong radial gradient of the ratio Mg2/Fe52 could reside in the 
different dependence on age for both indices. It is well known in fact that Mg2 changes with age 
more strongly than Fe52. This dependence is however already included in the synthesis models 
through the fitting functions which include a component that varies with the effective temperature 
of stars. Thus, the models should be able to predict the observational data. 

Maybe our models do not contain a sufficient variation in the mean age of the stellar populations 
along the radius. It is possible to obtain stronger radial gradients in age by decreasing the collapse 
time scale, for example TcoU = 1 Gyr, or by taking a shorter scale length X^. Both possibilities 
have been already included in our initial 500 models. But it should be noted that these types of 
models evolve very rapidly in the inner disk by swiftly consuming all the gas, even the molecular 
phase. The resulting radial distributions of molecular gas density show a hole in the center in 
clear disagreement with the observed exponential shape. It is also useful to remember that the 
uncertainties in the estimations of the molecular phase mass is large, due mostly to the dependence 
of its calibration on the metallicity of the region under study. We suspect that if this dependence 
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is included, the density of molecular gas in the inner disk could decrease and the models would 
change accordingly. 

Other possible solution is that the age dependence of Mg2 or Fe52 included in the fitting 
functions is inadequate. The fitting function for Mg2 from Borges et al. (1995) has a stronger 
dependence on age than the one from Worthey (1994), but the parameters range for the first library 
is narrower that for the second one. The latter may produce problems when young populations 
coexists with the older ones, as is the case in spiral disks. This is particularly so in the outer 
regions, precisely where models differ more from the observations. 

In our simple scenario the disk is created from the collapse of the gas of the protogalaxy. One 
might however think that other mechanism such as mergers are also able to produce variations in 
the ratio of Magnesium to Iron. Existing models for elliptical galaxies seem to indicate that in the 

merger scenario a ratio [Mg/Fe] > 0, and larger in the center, cannot be obtained (Thomas 2000). 
In fact the observed trends in [Mg/Fc] are the opposite to those derived from the semianalytical 
models simulating mergers of spiral to form ellipticals. 

Finally we must take into account the uncertainties associated with the difficult measurements 
of external disk regions of low surface brightness. It is clearly necessary to obtain a larger dataset 
of absorption spectral indices in spiral disks for galaxies with strong radial gradients of nebular 
abundances. 

4. Discussion: the possible evolutionary histories 

The impact of radial distributions of spectral indices on the model predictions can be explained 
through the effect of the star formation histories on the average stellar abundances. The star 
formation histories of our 19 models are different and this must imply different mean ages for the 
stellar populations along galactocentric radius. 

This may be clearly seen in Fig. 5 where we represent the star formation history for a inner 
(i? = 2 kpc) and a outer (R = 14 kpc) disk regions versus the oxygen abundance for two models. 
We have selected for comparison Models 129 and 450, which we consider clearly separated in the 
input parameter space. The first one is represented by the solid line while the dotted line shows 
results for the second one. The inner regions have triangles as symbols, and the outer squares. We 
see that Model 129 is above Model 450 for both regions. Although the final results are not very 
different, the evolutionary histories truly are, mostly in the inner regions. Therefore mean ages and 
abundances produced by these models must differ, and, as a result, the spectral indices must be 
different. 

We have computed these averaged ages and abundances for these two regions of the 19 models 
from Table 1. We represent our results in Fig. 6, panel a) for = 2 kpc and panel b) for R = lA 
kpc. 
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Fig. 5. Evolution of tlic star formation rate ^ vs the oxygen abundance for the inner [R = 2 
kpc, A symbol) and outer [R = 14 kpc, □ symbol) regions for Models 129 (solid lines) and 450 
(dotted lines). 
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Fig. 6. — The corresponding averaged age of the stellar populations vs the averaged oxygen abun- 
dance < [O/H] > for the 19 models in 2 radial regions: a) R=2 kpc and b) R= 14 kpc. The two 

vertical lines in panel a) join models with T = 4 (dot-long dashed line) and T = 4.5 (short-long 
dashed line). The square at the left , Model 450, is T = 5. Horizontal lines represent tq = 8 
(long-dashcd line), 12 (short-dashed) and 16 Gyr (dotted line), respectively. The top square is 
To = 4 Gyr. Symbols mark different = A^: A for = 4; □ for A^ = 8; • for A^ = 12 and ★ for 
Xd = 16 kpc, respectively. Panel b) Symbols have the same meaning than in Panel a). 
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In panel a), with a point for each model, we see an almost orthogonal diagram. The near- 
horizontal lines join models with a same collapse time scale, with different symbol following the 
scale length value (A for Aj^ = 4; □ for A/j = 8; • for A^ = 12 and -k for Xd = 15), while the vertical 
lines correspond to different efficiencies or T's. Spectral indices depend on both the abundances 
and the mean ages of the stars which produce them, and we see that Model 129, located at the top 
of the panel as a triangle, is older and more metal-rich that Model 450, the square located at the 
left of the same panel. The corresponding spectral indices can not be equal. The other models arc 
located in intermediate positions. The situation for outer regions ( i? = 14 kpc) is shown in panel 
b) where all models fall in a line. 

We can see, however, that for both panels the larger the mean abundances the older the stellar 
population. This means that the degeneracy problem of Mg2-Fe52 is not at work here because 
the result goes in the opposite direction. If the abundances turned out lower in a model than in 
other we would need a mean age older for the first in order to obtain similar spectral indices, while 
in fact our low-metal rich models turned out to be the youngest ones. Spectral indices, specially 
the Mg2 index, which has a strong dependence on mean age and has a better measurement, will 
probably discriminate among models better than the gas abundances. This effect is stronger in the 
outer regions of the disk, as we show in panel b), implying that spectral indices radial distributions 
may be similar in the inner regions but not in the external ones, where the gas distributions may 
diverge. This explains the effect shown in Fig. 4. 

The star formation histories have two effects on spiral disks: a radial gradient of abundances 
and a variable mean age for the stellar populations created in each radial region. Thus, even if the 
radial gradient of abundances obtained from different models are close enough to prevent adequate 
discrimination, a radial gradient of ages may still exist. Since spectral indices have a dependence 
on age and on abundances, a radial gradient in the mean ages of the stellar populations implies 
that radial distributions of indices for different models must differ more than those of abundances, 
thus helping in selecting the most adequate model or, at least, determining if some of these fall out 
of the confidence region of the plane defined by the data and their error bars. 

5. Conclusions 

We have discussed some chemical evolution models and their comparison with data. A first 
type of models, based purely on the multiphase chemical evolution models that yields direct nebular 
abundances were constrained by the radial distributions of both atomic and molecular gas. These 
observational data for the present time (Fig. 2), may be reproduced by only 19 among 500 initial 
models. Put differently, there arc some restricted sets of physically plausible input parameters that 
provide a superior fit to the abundance data. 

Next we computed the spectral index values predicted by the same models through the use 
of evolutionary synthesis models. These computations relay on a synthesis procedure capable of 
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establishing the abundance-index relation and which is implemented via the application of the first 
category of models above. This second category of models, which depends on the first one, may 
further constrain the chemical evolution of spiral disks by introducing measurements of their stellar 
content via the determination of spectral absorption indices, such as the Mg2 and Fe5270 Lick 
indices. Once again we show that the multiphase model reproduce the generic characteristics of 
both types of observed radial distributions. The 19 models initially selected on the basis of the 
nebular data do not reproduced equally well the radial distributions of spectral indices, as can be 
seen in Table 2: only 5 of them are now good enough to satisfy our criteria. 

We also show that statistical differences among models increase when spectral indices radial 
distributions are used with respect to those using only the present-day oxygen abundances. This 
effect appears because, even if the present abundances are similar for different models, the averages 
abundances and ages of the stellar populations responsible for the corresponding spectral indices 
are not. In fact, in a diagram abundance-age, our results go in the opposite direction than expected 
from the age-degeneracy problem. 

Therefore, by taking into account the results of Table 2, we demonstrate that spectral indices 
used as time constraints on the evolution of the galaxy are strong enough to eliminate a signifi- 
cant fraction of the models, or to help in the selection of possible scenarios, and we suggest that 
the measurement of spectral indices must be considered of importance in designing observational 
projects dealing with the structure of spiral disks. 
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